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Abstract 

An overview of the NASA ongoing efforts to 
develop an aircraft icing analysis capability is 
presented. Discussions are included of the over- 
all and long term objectives of the program as 
well as current capabilities and limitations of 
the various computer codes being developed. Des- 
criptions are given of codes being developed to 
analyze two and three dimensional trajectories of 
water droplets, airfoil ice accretion, aerodynamic 
performance degradation of components and complete 
aircraft configurations, electrothermal deicer, 
fluid freezing point depressant antideicer and 
electro-impulse deicer. The need for bench mark 
and verification data to support the code develop- 
ment is also discussed, and selected results of 
experimental programs are presented. 

Introducti on 

The aircraft icing problem has long been 
researched and studied. Reports exist in the 
literature which trace icing research activities 
back as early as the late 1920's. Since then 
research into the aircraft icing problem has been 
an almost ongoing effort at varying levels of 
intensity by a number of government and private 
research organizations both in the United States 
as well as in many other countries. A common 
thread which binds these various activities is 
that for the most part they have been experimental 
in nature. The majority of these programs have 
been aimed at studying the performance of various 
ice protection systems, although some limited 
attempts have been made to develop icing analysis 
capabi 1 ities. 

, The current NASA Aircraft Icing Research 
Program which was started in 1978 k 1 T seeks to 
take advantage of available computational fluid 
dynamics capabilities and develop a series of 
compatible computer codes which will address the 
fundamental icing problems. Once these initial 
codes are developed, additional codes will then 
be developed to broaden and enhance the capabili- 
ties. The computer codes will be thoroughly 
evaluated by comparison with appropriate verifi- 
cation data. The development of these computer 
codes and the acquisition of the required verifi- 
cation data bases are two major goals of the NASA 
aircraft icing research program. 

Once these various computer codes have been 
developed and adequately validated, they will 
provide the aircraft designers and certifiers 
with tools which should help to reduce the time 
and, therefore, the cost for the development and 
certification of aircraft with the required ice 
protection . 

Reference 2 presented the overall icing 
analysis plan which was developed and the status 


of the various computer codes at that time 
(1984). This report can be viewed as an update 
of that document. 
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Symb ol s 

aircraft drag and lift coefficients 

airfoil drag, lift, and moment 
coefficients 

airfoil chord 

icing cloud liquid water content, g/m^ 

icing cloud median volume droplet 
diameter, urn 

electrothermal deicer temperatures, °C 
free stream temperature, °C 
time, sec 

boundary layer velocity ratio 
free stream velocity, m/sec 
nondimensional coordinates 


a angle-of-attack 

Overall Objective^ 

The long range objective of the NASA aircraft 
icing analysis program is to develop a capability 
to predict the details of an aircraft icing 
encounter for both fixed and rotary wing vehicles. 
Figure 1 suggests four uses which could be made 
of such an analytical capability. Namely, the 
analysis could be used to predict aircraft per- 
formance and handling qualities changes due to 
icing on unprotected and de-iced components and 
to determine ice protection system performance. 

In addition, the analysis could be used to design 
ice protection systems for various applications 
and perhaps even be used to design components 
which are insensitive to icing. 

A majority of the current code development/ 
verification efforts are directed toward the two- 
dimensional airfoil problem. Once the capability 
to analyze the two-dimensional airfoil has been 
developed and the accuracy verified, it will serve 
as the basis for developing a fundamental analy- 
tical capability for treating wings, propellers, 
helicopter rotors, and complete aircraft con- 
figurations. 

Obviously a large number of computer codes 
are required in order to develop an overall 
aircraft icing analysis capability. Figure 2 
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shows the codes required and some (but by no 
means all) of the interfaces required to develop 
a unified methodology. 

The shaded boxes in Figure ? indicate those 
computer codes currently under development and 
verification. As the figure indicates, there are 
many computer codes which are not currently being 
developed. While available resources certainly 
limit the number of codes which can be developed 
concurrently, another, even more significant, 
limitation is a lack of knowledge of appropriate 
fundamental physics on which to base computer 
models. For example, the lack of understanding 
of basic structural and fracture properties of 
impact ice make it difficult to develop a code 
which models the pneumatic boot, one of the oldest 
of all ice protection systems. 

Thus, before computer codes can be developed 
to treat various aspects of the aircraft icing 
problem, a series of basic modeling experiments 
first must be conducted in order that the key 
physical elements be incorporated in any subse- 
quent mathematical models which are formulated. 
This work is currently underway in several areas 
such as determining the mechanical properties of 
impact ice.^'. 

Figure 2 also points out that an integral 
part of the icing analysis methodology under 
development is the aerodynamic analysis of both 
individual components and even complete aircraft 
configurations. Since these codes have histori- 
cally been developed for analysis of aerodynamic 
performance in nonicing environments, their use 
to study icing problems will require that sig- 
nificant modifications and improvements be made. 
Obviously, the more robust and accurate these 
codes are, the more robust and accurate the icing 
analysis methodology will be. 

The development of an aircraft icing analysis 
methodology requires that a large number of sup- 
porting benchmark (modeling) and verification 
experiments be run. The experiments require a 
wide range of ground and flight test facilities. 
Laboratory tests are required to do some of the 
most fundamental experiments while icing wind 
tunnel tests are needed to acquire data such as 
ice accretion, aero performance, and ice protec- 
tion system performance for selected components. 
Dry wind tunnel tests are needed to look at the 
detailed flowfield characteristics about airfoils 
with artificial leading edge ice accretions. 

Flight tests are required, both in simulated and 
natural icing environments, to determine overall 
aircraft performance and handling qualities 
changes . 

The following sections of the paper will 
discuss in detail the particular icing codes 
being developed as indicated in Figure 2. 

Pi scu ssion of Icing C o mputer Code s 

Water Droplet Traje ctory An al ysis 

A knowledge of the water droplet trajectories 
about any component provides the first indication 
of that component's sensitivity to icing. NASA 
is sponsoring the development of a family of drop- 
let trajectory codes which can handle a wide 


variety of geometries, from simple two-dimensional 
airfoils to complete aircraft configurations. 

Table 1 summarizes the various trajectory codes 
being developed and the geometries which can be 
handled. As the table suggests, the codes use 
various flowfield analysis methods. Two different 
techniques are used to pass flowfield information 
to the trajectory routine. The direct solution 
approach solves for the flowfield at each point 
along each trajectory and the grid generation 
approach calculates the velocity at the point of 
interest by interpolating on the two or three 
dimensional grid generated by the flowfield 
solver . 

There are advantages and disadvantages to 
both approaches. The direct solution approach 
can be more expensive computationally, especially 
as the number of trajectories computed increases. 
However, this approach frees the user from having 
to worry about generating a grid which adequately 
resolves the flowfield acceleration near the body 
surface. In addition, the direct solution approa- 
ches tend to have more difficulty computing the 
flowfield around a body such as an airfoil which 
has a leading edge ice accretion. Such calcula- 
tions are required as part of an ice accretion 
analysis as will be discussed later. The grid 
generation approach can handle airfoils with lead- 
ing edge ice accretions, but the development of 
an acceptable grid requires care and skill on the 
part of the analyst. 

All the trajectory codes indicated in Table 1 
use a predictor-corrector approach to integrate 
the droplet equations of motion. Different drop- 
let drag coefficient expressions are used by the 
various codes, but these differences are judged 
to not have a significant effect on the calcula- 
tion of droplet trajectories. 

Historically, water droplet trajectory codes 
have been evaluated by comparisons with available 
airfoil collection efficiency data obtained 
experimentally by NACA.^ 9 ' In general, the 
codes appear to be in reasonably good agreement 
with the experimental results which may themselves 
be subject to some appreciable errors. However, 
this database is limited in scope and in particu- 
lar, is confined to low speed studies of airfoil 
sections of interest to the aviation community in 
the 1940 to 1955 time period. Also, no experi- 
mental data had been available for inlet con- 
figurations of interest. 

An experimental program has been initiated 
to acquire a comprehensive collection efficiency 
database for a wide variety of modern airfoil and 
inlet geomerties which are of interest to the 
aviation community. The experimental approach, a 
variation on the dye tracer technique developed 
by NACA, is discussed in Reference 10. 

The major difference from the NACA technique 
is that a laser reflectance technique is used to 
determine the local collection efficiency distri- 
bution curves from the blotter paper strips which 
have been placed at various locations on the test 
model. The NACA approach was to dissolve out the 
dye from small samples of the blotter paper strips 
and then perform a colorimetric analysis of each 
resultant water-dye mixture. Figure 3 shows a 
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close-up photo of a quarter scale Boeing 737-300 
inlet model in the NASA Icing Research Tunnel 
(IRT) with a series of blotter strips attached at 
various circumferential locations. The laser 
set-up used to determine the collection efficiency 
distribution is shown in Figure 4. A Helium-Neon 
laser is used to provide the coherent light source. 
The amount of light reflected from the blotter can 
be related to the local amount of dye on the blot- 
ter. The dye was dissolved in the spray water 
prior to being feed through the spray nozzles. 

This dye concentration is directly related to 
local collection efficiency. This automated data 
reduction procedure allows a much larger quantity 
of data to be reduced relative to the labor inten- 
sive method which had to be employed by NACA 
researchers . 

Figure 5 shows measured collection efficiency 
"data for the o NACA 65?-0l5 airfoil for two angles 
of attack (0° and 8°). Also shown are predicted 
collection efficiency curves. The results are in 
close agreement for the 0° angle-of-attack case 
while some differences are noted for the 8° case. 

Initially, drop impingement tests will be 
conducted in the NASA Icing Research Tunnel (IRT) 
using large number of airfoil and inlet models. 

Once these tests are completed, follow-on tests 
will be conducted in a high speed icing wind tun- 
nel to quantify the effects of flowfield com- 
pressibility on droplet trajectories. Previous 
analytical studies conducted jointly by NASA and 
the British RAE^l) have indicated that com- 
pressibility becomes important only for the small- 
est dropsizes. However, this analytical result 
must be verified experimentally. 

Some current applications of the trajectory 
analysis codes will be mentioned. Figure 6 shows 
selected results of three-dimensional trajectory 
analysis(°) for the nose region of the NASA 
icing research aircraft, a deHavilland Twin Otter. 
The purpose of this study was to look at the 
effects of the flowfield on droplet trajectories 
and therefore, on the measurements of icing cloud 
liquid water content as measured by various in- 
struments located on the aircraft. The results 
are shown in the form of concentration factor as 
a function of droplet diameter. Concentration 
factor indicates the number of droplets per unit 
volume at the measurement station ratioed to the 
number of droplets of that size per unit volume 
in the free stream cloud. As the figure indi- 
cates, for the selected location, the instrument 
sees more water droplets per unit volume than 
actually exist in the cloud especially, for a 
range of droplet sizes between about 20 and 200 urn 
and thus would indicate a liquid water content 
somewhat higher than the freestream value. This 
analysis is being used to help in interpreting 
icing cloud data currently being acquired and to 
assist in future placement of instruments. 

An indepth study of airfoil collection 
efficiency characteristics has been investigated 
using the two-dimensional trajectory code of 
Reference 4. Some 30 different airfoil geometries 
were analyzed for a wide variety of velocity and 
angle-of-attack conditions for a range of water 
droplet sizes. The results were used to form a 
computational database to investigate the effects 
of airfoil design parameters such as leading edge 


radius, camber, and maximum thickness on collec- 
tion efficiency .U2) j n addition, correlation 
expressions have been developed to allow quick 
estimates to be made of airfoil collection effic- 
iency characteristics for geometries other than 
those analyzed. 

The three-dimensional trajectory analysis 
code discussed in References 6 and 7 is being 
applied to inlet configurations which have flow- 
fields which are highly three-dimensional in * 
nature. The initial results ( ^ ) are highly 
encouraging and indicate complex inlet geometries 
can be analyzed for icing sensitivity. 

Airfoil I ce Accretion 

The NASA airfoil ice accretion analysis code 
which is called LEWICE is an extension of the 
work of Lozowski et al.”^' and Ackley and 
Templeton. U5) The LEWICE code predicts the 
ice growth rate distribution (and thus ice shape) 
around the leading edge of the airfoil by locally 
solving the quasi-steady energy balance equation 
first proposed by Messi nger . ' The energy 

balance accounts for the governing heat and mass 
transfer processes thought to occur during the 
icing process. These terms include convection to 
the free stream, latent heat release due to freez- 
ing, sublimation (or evaporation), and aerodynamic 
heating. A more detailed discussion of the LEWICE 
code is given in Reference 17. 

The four major elements of the LEWICE code 
are (1) a potential flow panel code to calculate 
the airfoil flowfield, (21 a two-dimensional water 
droplet trajectory code,^ 5 ) (3) an energy bal- 
ance routine, and (4) integral boundary layer 
routine to predict heat transfer distributions. 

The LEWICE code is constructed so that the flow- 
field and droplet trajectory calculations are 
repeated at user specified time intervals as the 
ice accretion grows on the airfoil leading edge. 
This updating process attempts to account for the 
effect of the leading edge ice accretion on the 
airfoil flowfield and thus, on the changes in 
airfoil collection efficiency characteristics. 

The current capabilities of the LEWICE code 
are summarized in Figure 7 where ice accretion 
shapes measured on a NACA0012 airfoil in the IRT 
are compared with LEWICE predictions. Two sets 
of comparisons are shown for 5 min rime and glaze 
ice accretions. For each calculation, the flow- 
field and droplet trajectory analyses were updated 
every 1 min. The agreement is judged to be good 
for the rime ice case, but the LEWICE code over 
predicts the glaze ice accretion, although the 
general shape appears to be good. 

The integral boundary layer code has a simple 
treatment of the effects of surface roughness on 
boundary layer transition and heat transfer levels. 
However, this model is formulated in terms of 
equivalent sand grain roughness and the relation- 
ship to actual ice surface roughness is still a 
research topic, so currently, "educated guesses" 
must be input to the model. Other models for 
wall roughness currently are being considered for 
incorporation in the LEWICE code. 

Recent high-speed photo-micrographic movies 
of the icing process have suggested that 
some of the fundamental assumptions made in the 
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ice accretion model may be in error. Efforts 
will shortly begin to formulate either changes to 
the existing model or an alternate model in order 
to see if ice shape predictions can be improved. 

Other research begun recently in order to 
support the ice accretion modeling efforts 
include the measurement of convective heat 
transfer coefficients from a smooth and rough 
surface airfoil both in low and high turbulence 
environments and application of a Navier-Stokes 
solver to predicting heat transfer coefficient 
distributions around bodies with ice accretions. 

Aerodyna mic Pe rformance Degradation Analysis 

Accurate predictions of aerodynamic perform- 
ance degradation of an aircraft due to icing is 
one of the desired end products of the icing 
analysis methodology. Currently, the development 
and verification of computer codes to predict 
airfoil performance degradation in icing is the 
primary research effort, but some more empiri- 
cally based approaches are being developed for 
predicting degradations of propellers, helicopter 
rotors, and complete aircraft. The analysis for 
the airfoil will provide the basis for developing 
more exact analyses for components such as wings, 
tails, propellers, rotors, and eventually for 
complete aircraft configurations. 

Currently, two advanced analysis codes are 
being evaluated for use in predicting iced airfoil 
performance - a Navier-Stokes solver and an inter- 
active boundary layer code. Code predictions are 
being compared to data being gathered in the IRT 
and in the Ohio State University 0.9 by 1.5 m 
low-speed wind tunnel. These codes were chosen 
since they have both demonstrated the ability to 
treat flows which have regions of separation- 
reattachment, a phenomenon which often occurs 
with ice shapes. Previously, more conventional 
a i rf o i 1 analysis codes were considered, but the 
codes all were incapable of handling the iced 
airfoil problem without major revisions being 
i ncorporated. A summary of these investigations 
is given in Reference 2. 

1 he Navier-Stokes code being evaluated is 
the ARC2D code developed the at NASA Ames Research 
Center. ^ iy ' The code has been used to analyze 
the flowfields about a large number of airfoil 
geometries including airfoils with extended spoil- 
ers which create a large separated zone. The code 
can be run in either of two modes - the full 
Navier-Stokes mode or the so-called thin layer 
Navier-Stokes mode which removes the viscous terms 
in the axial direction. To date, most calcula- 
tions have been done with the thin layer Navier- 
Stokes mode, although, limited comparisons of the 
predictions made with both modes failed to reveal 
any significant differences. The ARC2D code uses 
the Baldwin-Lomax turbulence model, and the so- 
called GRAPE grid generation code developed by 
Sorenson, also of NASA Ames.^O) 

Initial ARC2D calculations were made for a 
NACA632-A415 airfoil with a simulated leading 
edge glaze ice shape. The ARC2D results 
are compared in figure 8 to measurements made in 
the IRT with a general aviation wing section with 
a NACA632-A415 airfoil section with a wooden 
leading ice shape affixed. The predicted results 
(assuming no laminar boundary layer growth) are 


judged to be in good agreement with the IRT data 
with the exception of the over prediction of 
Ci . Similar over predictions of Cl with 
ARC2D have been found by other researchers for 
clean airfoil configurations. The agreement 
in drag polars is thought to be exceptionally 
good. The velocity vector plot shown in the fig- 
ure indicates that the ARC2D code predicted siz- 
able zones of separated flow on both upper and 
lower surfaces, but reattachment was aslo pre- 
dicted for both surfaces. t 

While these initial results suggested that 
the ARC2D code could predict the performance of 
airfoils with leading edge ice shapes, the data 
acquired in the IRT was judged to be not detailed 
enough to serve as a val idation/verif ication data 
base. It was determined that such a database 
should include, in addition to the airfoil per- 
formance data (Ci, Cd, C m ) , much more detailed 
surface pressure distributions, velocity and tur- 
bulence profiles, especially, in the separation- 
reattachment zones, and flow visualization data. 

In order to acquire this database, an extensive 
experimental program was initiated in the Ohio 
State University low speed wind tunnel., Two 
summary papers of that effort are given in 
References 23 to 24. 

The airfoil geometry chosen for the program 
was a 0.53 m (21 in.) chord NACA0012 airfoil. A 
removable leading edge was fabricated which had 
an ice shape which approximated an ice accretion 
measured in the IRT. The ice shape and its rela- 
tion to the IRT ice shape are shown in Figure 9. 
The ice shape geometry can be described by a con- 
stant radius section, two rounded edges of speci- 
fied radius, and two straight line sections. 

Thus, an exact specification could be used to 
input the geometry into any flowfield analysis 
code. The wind tunnel model was instrumented 
with a surface static pressure every 1 percent 
chord in the leading edge region. This large 
number of pressures gave excellent definition of 
the separation-reattachment zones on the surfaces. 

Figure 10 presents the ARC2D predictions 
compared to the OSU data. Again, the agreement 
is judged to be quite good for the lower angles- 
of-attack. However, it was observed that the 
ARC2D code did not converge to a steady state 
solution for angles of attack greater than about 
6°. Currently, flow visualization experiments 
are being conducted to determine if this is purely 
a computational result, or if the experimental 
flowfield is not steady at these higher angles- 
of-attack . 

To date, velocity and axial turbulence 
intensity profiles have been measured every 
2 percent chord for the upper surface separation- 
reattachment zone. Figure 11 shows two selected 
profile comparisons - one just downstream of the 
measured point of separation and the other just 
upstream of the experimental reattachment point. 
The ARC2D predictions are in reasonable agreement 
with the data; however, the predicted reattachment 
occurred more upstream than was observed in the 
experiment. Investigations are underway to 
improve the code predictions. Areas of investi- 
gation include the sensitivity to ice shape defi- 
nition and grid generation as well as to the 
turbulence model employed. 
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A similar evaluation of the capabilities of 
the interactive boundary layer code was begun 
more recently. This code was developed by Cebeci 
and associates under grant to the National Science 
Foundation^ ) to predict clean airfoil behavior 
including the higher angles-of-attack and the 
onset of stall. The code has two major modules - 
a potential flow panel routine to calculate the 
inviscid region and an inverse finite difference 
boundary layer routine which is capable of treat- 
ing regions of separation-reattachment. The two 
routines are coupled together to solve further 
airfoil flowfield by iterating back and forth 
until a converged solution is acquired. 

Initial studies have been completed for lower 
angle-of-attack results for the "iced" NACA0012 
airfoil already discussed with comparisons made 
'of airfoil performance (i.e., Ci, C^). The pre- 
dictions were in very close agreement with the 
experimental data. Some oscillations in the final 
solutions obtained suggested some additional 
improvements are required in the numerics. These 
improvements are currently being made and the 
more detailed velocity profile comparisons are 
also being made. Also, comparisons are being 
made for the higher angle-of-attack conditions. 

Currently, it is judged that both the Navier- 
Stokes and interactive boundary layer methods 
show great promise for being able to predict iced 
airfoil performance. It is felt that the two 
methods are more complementary than competitive 
as both can be used for icing analysis. The 
Navier-Stokes approach, while computationally 
more expensive (approximately 10 to 20 min per 
solution on the Cray XMP) is thought to be capable 
of treating larger, more extreme ice shapes. 
Navier-Stokes results, possibly, can be used to 
help in the more approximate modeling in the 
interactive boundary layer code. The interactive 
boundary layer code is much faster (approximately 
4 to 30 sec per solution on the Cray XMP) and 
would be very attractive for those problems which 
require many aero performance calculations to be 
made (e.q. the helicopter rotor in forward 
flight). 

Besides those areas of improvements currently 
being investigated for both codes which have been 
discussed, each will require a model to be devel- 
oped to handle the extreme levels of surface 
roughness which exist for any ice shape. At pre- 
sent, neither code can handle anything other than 
smooth surfaces. Several approaches for treating 
surface roughness are being evaluated. Regardless 
of the model employed, some means must be found 
for quantifying the nature of the ice shape 
roughness. Currently, the equivalent sand grain 
approach is employed to treat rough surfaces, but 
it is not clear that that approach would work for 
ice accretion roughness. 

As already indicated, more approximate meth- 
ods have been developed for predicting the aero 
performance degradation of propellers, helicopter 
rotors (hover and forward flight), and complete 
aircraft. References 26 to 30 review these 
efforts. These various efforts all use the same 
approach - a set of correlation equations is used 
to predict airfoil performance degradation due to 
icing and, then an appropriate performance code 


is used to determine the total component or air- 
craft performance degradation. 

The correlation equations are relations 
developed from an airfoil icing data base which 
relate the change in lift, drag, and pitching 
moment to known aerodynamic and icing variables. 

The first such correlations were developed by 
Gray of NACA.( 31 ' More recently, Bragg devel- 
oped a rime ice correlation.^ ' Both of these 
correlations used data gathered in the NASA FRT . 
Flemming, (29) acquired a large data base in the 
Canadian NRC high speed icing wind tunnel for a 
series of reduced scale rotor airfoil sections 
which he used to develop a series of correlation 
equations. 

Unfortunately, it is not possible to a 
priori select one set of correlation equations 
with assurance that the predictions will yield 
accurate results. Studies made with the above 
mentioned sets of correlations have determined 
that the errors can be as large as 100 to 
200 percent for selected conditions. Some of 
this error is undoubtedly due to the rather limi- 
ted airfoil icing data base which exists. Also, 
questions can be asked about how the icing results 
from one facility can relate to those from another 
facility which were calibrated using different 
techniques. Thus, extreme care must be exercized 
when using these correlations . Miller^! gives 
an excellent discussion of typical effects of 
using the different correlation equations to 
predict propeller performance. 

In spite of the above concerns, some reason- 
able results have been acquired. Figure 12 shows 
propeller performance predicted in rime icing 
conditions when compared to natural icing data 
gathered by NACA researchers. The propeller effi- 
ciency variation as a function of advance ratio 
measured in flight icing conditions was closely 
duplicated by the analytical methodology for this 
particular set of tests. As Reference 26 indi- 
cates, the agreement was not as good for some of 
the other rime icing encounters. 

Korkan, et al., used a similar analytical 
approach to analyze the helicopter rotor, first 
in hover(2°) and then in forward flight. (2- ' 

The particular rotor configuration analyzed was 
the front rotor of the Boeing Chinook CH47 . 
Unfortunately, no experimental data existed for 
comparison with the predictions. All that could 
be said was the results were reasonable when com- 
pared with previously reported torque rises mea- 
sured in icing for other heliocpter configurations. 

Flemming used the previously discussed air- 
foil performance in-icing correlations along with 
appropriate rotor performance codes to predict 
rotor degradation for the S-76, UH1H, and UH-60A 
helicopters. The predicted results compared to 
available experimental data are shown in Fig- 
ure 13. Generally the agreement was better for 
the hover comparisons than for the forward flight. 
As Figure 13 indicates, the forward flight pre- 
dictions were generally lower than the experimen- 
tal measurements. At present,^ it is not known 
why this is true, although, as Flemming points 
out, (29) i t is not clear how accurate are the 
experimental measurements. These results point 
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out the need for dedicated icing flight programs 
for acquiring high quality verification data. 

Figure 14 shows selected results from the 
prediction of overall aircraft performance degrad- 
ation in icing using the NASA code^ u J compared to 
data gathered as part of the NASA icing research 
flight program^) using the Twin Otter. The 
airfoil correlation developed by Bragg^ 4 ' was 
used to determine airfoil performance degradation. 
Even though this correlation was developed to 
handle rime icing conditions, the figure indicates 
it appeared to provide reasonable inputs for the 
mixed icing conditions of the icing flight shown. 
The code, as currently configured, can predict 
aero performance degradations only for lifting 
components (e.g. wings, tails, struts). The con- 
tributions of other components such as wheels, 
flap' hinges, antennas, etc., must be estimated 
and input. For the predictions shown in Figure 14, 
this contribution was estimated by forcing the 
code to agree with the drag increase measured for 
the completely iced aircraft for the lowest lift 
coefficient data point (filled symbol in the fig- 
ure). This amount of drag increase was not var ied 
but kept constant for all the other code predic- 
tions. Note the reasonable agreement between 
theory and experiment not only for the fully iced, 
but also for the wing and wing plus tail deiced 
data. Additional comparisons of the code predic- 
tions using the various airfoil correlations 
available with additional natural icing flight 
data are planned. Also, simple methods for esti- 
mating the drag increase due to nonlifting com- 
ponents are being considered. 

Ice Protectio n S y stems 

[\|A3A "has^sponsored the development of a 
series of computer codes which model the electro- 
thermal deicer, a system currently being used for 
many helicopter rotor and fixed wing aircraft 
applications. Table 2 gives the important charac- 
teristics of the six different codes developed to 
date and how they differ from one another. As 
the table suggests, the codes vary. greatly in 
complexity and therefore, in the size of computer 
required. Codes 1 and 2 which model the electro- 
thermal deicer on a one-dimensional, time varying 
basis can be run on a personal computer. Code 6 
which analyzes the complete two-dimensional geo- 
metry, including the variable thickness ice layer 
and the many layers of the airfoil geometry 
requires a supercomputer (i.e. a Cray IS or XMP) 
to achieve reasonable run times. A unique feature 
of all codes is that they treat the melting of 
the ice through a phase change routine. Codes 2 
through 6 use the weak enthalpy method^- ' 
while code 1 uses an approximate method. 

Currently, the code predictions are being 
compared with electrothermal deicer data from 
icing wind tunnel and helicopter natural icing 
flight tests. These comparisons are aimed at 
providing a better understanding of how accu- 
rately the various codes can predict electro- 
thermal deicer performance for various levels of 
modeling complexity. 

Figure 15 shows the test model used in a 
recent IRT test. The model was a section of a 
Bell Flelicopter Textron UH1H rotor (0.53 m (21 in. 
chord)) with a spanwise electrothermal deicer 
configuration installed. The deicer, provided by 


the BF Goodrich Company, is typical of a flight 
design. The model was instrumented with a large 
number of thermocouples to measure the transient 
temperature response at various locations in the 
model. The test model was run over a wide variety 
of dry, wet (above freezing), and icing condi- 
tions, and a large data base was acquired. 
Reference 36 gives a description of the test 
program. 

Figure 16, taken from Reference 38, shows * 
one selected result of the comparison of code 2 
predictions with experimental measurements. As 
the figure shows, thermocouples were located at 
(1) the abrasion shield surface, (2) the heater, 
and (3) the internal spar of the model. For this 
comparison, the model predictions were judged to 
be in close agreement with the experimental data 
for all three cycles of electrothermal deicer 
operation. The model predicted a melting of the 
ice for the second and third heater cycles, and 
the experimental results showed the same occur- 
rence. While the agreement was not as good for 
all cases as the one shown here, the one- 
dimensional code predictions were found to be 
generally in good agreement with the data. 
Currently, efforts are underway to make similar 
comparisons with the other codes to see if the 
increased complexity in the numerical modeling 
can give better agreement with the experimental 
results. 

Ice shedding is a key aspect of electro- 
thermal deicer performance. The one-dimensional 
model of Reference 34 attempts to model this 
phenomenon by assuming that shedding occurs when 
the thickness of melted ice reaches a user speci- 
fied amount. A better treatment of shedding may 
be required. 

Another ice protection system, the fluid 
freezing point depressant system, has been 
modeled. A simple engineering procedure for pre- 
dicting the minimum anti-icing flowrates required 
has been developed.^ 9 ' The code predictions 
compared with flowrates measured in the IRT for 
a full-scale general aviation wing section model 
equipped with this system are shown in Fig- 
ure 17. In general, the method is judged to give 
reasonable agreement with the experimental data. 
As Reference 39 points out, some of the error can 
be attributed to the experimental values measured 
for minim um anti-icing flowrates. This deter- 
mination based on visual observations was found 
to sometimes vary significantly for repeat test 
conditions and was also sensitive to the inter- 
pretations of different observers as to when the 
minimum flow rate was achieved. Given these 
inherent uncertainties in the experimental 
approach, the predictive method is felt to be 
acceptably accurate. 

The third ice protection system being 
modeled is the electromagnetic impulse deicer 
(EIDI). EIDI is a concept which is receiving 
much attention, both in the United States and 
Europe. It shows great promise for being a low 
power consumption, highly efficient deicing sys- 
tem. NASA has sponsored an intensive research 
effort, the results of which are summarized in 
Reference 40, to acquire an experimental data 
base to aid in the design of the system for 
various aircraft applications. 
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Work is also underway to develop the required 
structural and electrodynamic codes to model the 
key features of the EIDI system. The current 
status of the structural dynamic modeling capa- 
bilities is discussed in Reference 41. The simple 
configuration which is being modeled is shown in 
Figure 18. The half cylinder is meant to model a 
wing leading edge and has a single 2.54 cm radius 
coil located 18 cm (0.07 in.) beneath the surface. 
Finite element model predictions for surface 
acceleration at two locations and the circumferen- 
tial and longitudinal strains at the coil location 
are shown in the figure. Agreement in all cases 
is judged to be at least reasonable and surpris- 
ingly good in some instances (e.g. peak circum- 
ferential strain). 

A key input into the structural dynamic 
analysis is the force produced by the coil. Both 
ther temporal and spatial characteristics of this 
force must be known. Until recently, no analyti- 
cal solutions had been obtained for this force 
function, and thus experimental data had to be 
used. The results shown in Figure 18 were 
acquired using a simple force model. However, 
recently, a transmission line model obtained from 
the governing field theory equations has been 
achieved which appears to yield good predictions 
of the force for the case of a coil with its axis 
perpendicular to a flat plate. '^) Figure 19 
shows a comparison of the force versus time pre- 
dicted by this transmission line model with the 
force versus time calculated from the measured 
induction fields for a 3.08 cm (2 in.) diameter 

0.478 cm (0.188 in.) thick coil rigidly mounted 
next to a fixed 0.0813 cm (0.032 in.) thick 2024 
T3 Aluminum disc with a 12.7 cm (5 in.) diameter. 
It is hoped that the use of this force model will 
improve theoretical predictions such as those 
shown in Figure 18. 

While the modeling of the EIDI system shows 
great promise, the effort is still in its infancy. 
A better understanding of the ice adhesion prop- 
erties and removal mechanisms must be gained 
before a complete modeling of the system can be 
attempted. Also, the EIDI system modeling for 
most, applications will probably require three- 
dimensional structural analysis codes. 

As already indicated, NASA has sponsored 
efforts to measure some of the key physical 
properties of ice^) which are required as 
inputs into any surface deflection type ice 
protection system model. The results are being 
used along with a finite element analysis approach 
to begin a quasi steady modeling of a generic 
pneumatic boot ice protection system. Companion 
experiments are being planned to measure the ice 
fracture process of the generic boot configuration 
for comparison with modeling predictions. 

Concluding Remarks 

This paper has presented a current status 
review of the progress being made toward devel- 
oping and validating the computer codes required 
to assemble an aircraft icing analysis method- 
ology. While it is felt that much progress has 
been made since the initial review presented in 
Reference 2, much work remains to be done. While 
many of the codes which are needed to form a core 



capability have been developed, many validation/ 
verification experiments are still required. This 
is especially difficult (and expensive) when it 
comes to determining the performance in icing of 
complex configurations such as helicopter rotors 
and complete aircraft configurations. Analysis 
of the steady state and dynamic behavior of those 
complex configurations in icing may require sig- 
nificant levels of empiricism to be incorporated 
i n the modeling. 

t 

Nevertheless, it is felt that once developed 
and validated, this icing analysis capability will 
be of great benefit to the industry in reducing 
the fixed and rotary wing development and certi- 
fication costs for future aircraft configurations 
which are required to safely fly into forecast 
icing conditions. 
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TABLE 1. - DROPLET TRAJECTORY ANALYSIS CODE CHARACTERISTICS 


Geometries 

Flowfield analysis code used 

Comments 

Reference 

1. Single Element Airfoils 

Theodorsen Transformation 

Incompressible Flowfield, 
Direct Solution 

4 

2. Single Element Airfoils 

Navier-Stokes 

Compressible Flowfield, 
Grid Generation 

5 

3. Multiple Element Airfoils 

Panel Method 

Incompressible Flowfield, 
Direct Solution 

5 

4. 2D, Axisymmetric Inlets 

Panel Method 

Incompressible Flowfield, 
Direct Solution, 

2D Trajectories 

5 

5. 3D Inlets 

Transonic, Full Potential 

Compressible Flowfield 
Grid Generation 

6 

6. Complex Configurations 
(incl. Complete Aircraft) 

Panel Method 

Incompressible Flowfield, 
Direct Solution 

7,8 


TABLE 2. - ELECTROTHERMAL DEICER CODES 


Geometry 

Ice thickness 

Number of heaters 

Phase change 

Reference 

1. One Dimensional 

Constant 

One 

Approximate 

33 

2. One Dimensional 

3. Two Dimensional 

Constant 

Constant 

One 

One 

Weak Enthalpy 
Weak Enthalpy 

34 

35 

(Recti linear) 




4. Two Dimensional 

Constant 

Vari abl e 

Weak Enthalpy 


(Recti 1 inear ) 




5. Two Dimensional 

Variable 

One 

Weak Enthalpy 

37 

(Rectilinear) 

6. Two Dimensional 

Vari abl e 

Multiple 

Weak Enthalpy 

36 

(Exact ) 
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Figure 2. - Aircraft icing analysis methodology. 
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Figure 5. - Preliminary collection efficiency results for NACA 
652~015 airfoil. V*, = 82 m/sec ; MVD = 15 |im. 



Figure 6. - Concentration factor calculations for NASA icing 
research aircraft, = 67 m/sec. 
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' Figure 7. - LEWICE predictions of airfoil ice accretions compared to IRT data. NACA 0012, 0.53 m 
chord, 5 minute ice accretion, 0° angle of attack. 
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Figure 8. - Navier-Stokes analysis 
of NACA 632 _ A415 airfoil with 
glaze ice (ref. 21). 
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Figure 9. - Comparison of measured and simulated ice shapes 
NACA 0012 airfoil, 0.53 m <21 in. ) chord. 






DIVIDING 

STREAMLINE 

(EXPERIMENTAL)- 



Figure 11. - Comparison of measured and predicted velocity profiles using ARC 2D 
code. 
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Figure 12. - Comparison of measured and predicted 
propeller performance in rime icing conditions 
(ref. 26). 
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(b) Forward flight. 

Figure 13. - Comparisons of measured and predicted 
helicopter performance in natural icing (ref. 29). 
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Figure 14. - Comparison of predicted and measured aircraft per- 
formance in natural icing conditions (ref. 30). 



Figure 15. - Electrothermal de-icer experiment in IRT. 
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Figure 17. - Comparison of predicted and measured 
minimum anti-icing flow rates for NACA 2412 
(mod) airfoil with fluid freezing point ice protec- 
tion system (ref. 39). 
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Figure 18. - Comparison of measured and predicted accelerations and strains for idealized deicer model 
(ref. 41). 
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Figure 19, - Coil-metal target separation force. 
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